﻿/********************************************************
 *  ██████╗  ██████╗████████╗██╗
 * ██╔════╝ ██╔════╝╚══██╔══╝██║
 * ██║  ███╗██║        ██║   ██║
 * ██║   ██║██║        ██║   ██║
 * ╚██████╔╝╚██████╗   ██║   ███████╗
 *  ╚═════╝  ╚═════╝   ╚═╝   ╚══════╝
 * Geophysical Computational Tools & Library (GCTL)
 *
 * Copyright (c) 2023  Yi Zhang (yizhang-geo@zju.edu.cn)
 *
 * GCTL is distributed under a dual licensing scheme. You can redistribute 
 * it and/or modify it under the terms of the GNU Lesser General Public 
 * License as published by the Free Software Foundation, either version 2 
 * of the License, or (at your option) any later version. You should have 
 * received a copy of the GNU Lesser General Public License along with this 
 * program. If not, see <http://www.gnu.org/licenses/>.
 * 
 * If the terms and conditions of the LGPL v.2. would prevent you from using 
 * the GCTL, please consider the option to obtain a commercial license for a 
 * fee. These licenses are offered by the GCTL's original author. As a rule, 
 * licenses are provided "as-is", unlimited in time for a one time fee. Please 
 * send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget 
 * to include some description of your company and the realm of its activities. 
 * Also add information on how to contact you by electronic and paper mail.
 ******************************************************/

#ifndef _GCTL_TRIANGLE_H
#define _GCTL_TRIANGLE_H

#include "vertex.h"
#include "entity.h"

namespace gctl
{
	// Declaration of the basic triangle type
    template <typename A> struct type_triangle;
	// triangle type of attribute type of void
    typedef type_triangle<void> triangle;

	/**
	 * @brief      Structure of a triangle.
	 * 
	 * 三维直角坐标中的三角形结构体。
	 */
    template <typename A>
	struct type_triangle : public entity<vertex3dc, 3, A>
	{
		/**
		 * 构造函数
		 */
		type_triangle();
		/**
		 * @brief      Set object from parameters
		 *
		 * @warning    This function will locate memories to store vertice
		 *
		 * @param[in]  p0     The first point
		 * @param[in]  p1     The second point
		 * @param[in]  p2     The third point
		 * @param[in]  index  The element index
		 */
		type_triangle(const point3dc &p0, const point3dc &p1, const point3dc &p2, int index = 0);
		/**
		 * @brief      构造函数重载
		 *
		 * @param[in]  index     索引值
		 * @param      vert0  顶点0的指针
		 * @param      vert1  顶点1的指针
		 * @param      vert2  顶点2的指针
		 */
		type_triangle(vertex3dc &vert0, vertex3dc &vert1, vertex3dc &vert2, int index = 0);
		/**
		 * @brief      析构函数
		 */
		virtual ~type_triangle(){}
		/**
		 * @brief      Set object from parameters
		 *
		 * @warning    This function will locate memories to store vertice
		 *
		 * @param[in]  p0     The first point
		 * @param[in]  p1     The second point
		 * @param[in]  p2     The third point
		 * @param[in]  index  The element index
		 */
		void set(const point3dc &p0, const point3dc &p1, const point3dc &p2, int index = 0);
		/**
		 * @brief      设置三角形的索引与顶点列表
		 *
		 * @param[in]  index     索引值
		 * @param      vert0  顶点0的指针
		 * @param      vert1  顶点1的指针
		 * @param      vert2  顶点2的指针
		 */
		void set(vertex3dc &vert0, vertex3dc &vert1, vertex3dc &vert2, int index = 0);
        /**
         * @brief      Calculate the triangle's area
         *
         * @return     area
         */
        double area();
		/**
		 * @brief      Calculate the local coordinates of the triangle
		 * 
		 *      (hgt_x, hgt_y)
		 *       / \
		 *      / | \
		 *     /  |  \
		 *    /   |   \
		 *   /---------\
		 *  (0,0)  (len, 0)
		 * 
		 * @param len    bottom length the triangle which equals to the longest edge length
		 * @param hgt_x  local x coordinate of the height of the triangle
		 * @param hgt_y  local y coordinate of the height of the triangle
		 * @param c_id   local index of the vertex ar the origin
		 */
		void local_coordinates(double &len, double &hgt_x, double &hgt_y, size_t &c_id);
	};

	template <typename A>
	type_triangle<A>::type_triangle() : entity<vertex3dc, 3, A>::entity(){}

	template <typename A>
    type_triangle<A>::type_triangle(const point3dc &p0, const point3dc &p1, const point3dc &p2, int index) : type_triangle()
    {
        set(p0, p1, p2, index);
    }

	template <typename A>
	type_triangle<A>::type_triangle(vertex3dc &vert0, vertex3dc &vert1, 
		vertex3dc &vert2, int index) : type_triangle()
	{
		set(vert0, vert1, vert2, index);
	}

	template <typename A>
    void type_triangle<A>::set(const point3dc &p0, const point3dc &p1, const point3dc &p2, int index)
    {
        if (index < 0)
        {
            throw out_of_range("Invalid index number, From type_triangle::set(...)");
        }

        for (int i = 0; i < 3; ++i)
        {
            this->vert[i] = new vertex3dc;
        }
        this->self_host = true;

        this->id = index;
        this->vert[0]->set(p0, 0);
        this->vert[1]->set(p1, 1);
        this->vert[2]->set(p2, 2);
        return;
    }

	template <typename A>
	void type_triangle<A>::set(vertex3dc &vert0, vertex3dc &vert1, 
		vertex3dc &vert2, int index)
	{
	    if (index < 0)
	    {
	        throw out_of_range("Invalid index number, From type_triangle::set(...)");
	    }

		this->id = index;
		this->vert[0] = &vert0;
		this->vert[1] = &vert1;
		this->vert[2] = &vert2;
		return;
	}

    template <typename A>
    double type_triangle<A>::area()
    {
        double edge_len[3] = {0.0, 0.0, 0.0};
        edge_len[0] = distance(*this->vert[0], *this->vert[1]);
        edge_len[1] = distance(*this->vert[1], *this->vert[2]);
        edge_len[2] = distance(*this->vert[2], *this->vert[0]);

        // calculate element's area
        double p = 0.5*(edge_len[0] + edge_len[1] + edge_len[2]);
        return sqrt(p*(p-edge_len[0])*(p-edge_len[1])*(p-edge_len[2]));
    }

	template <typename A>
    void type_triangle<A>::local_coordinates(double &len, double &hgt_x, double &hgt_y, size_t &c_id)
	{
		double edge_len[3] = {0.0, 0.0, 0.0};
        edge_len[0] = distance(*this->vert[0], *this->vert[1]);
        edge_len[1] = distance(*this->vert[1], *this->vert[2]);
        edge_len[2] = distance(*this->vert[2], *this->vert[0]);

		if (edge_len[1] >= edge_len[0]) c_id = 1;
		else c_id = 0;

		if (edge_len[2] > edge_len[c_id]) c_id = 2;

        // calculate element's area
        double p = 0.5*(edge_len[0] + edge_len[1] + edge_len[2]);
        double area = sqrt(p*(p-edge_len[0])*(p-edge_len[1])*(p-edge_len[2]));

		len = edge_len[c_id];
		hgt_y = 2.0*area/len;
		hgt_x = sqrt(edge_len[(c_id+2)%3]*edge_len[(c_id+2)%3] - hgt_y*hgt_y);
		return;
	}
}

#endif // _GCTL_TRIANGLE_H